syms x y
eqns = [x^2 + y^2 == 5, x - y == 1];
vars = [x y];
[x, y] = solve(eqns, vars)